setwd("")

rm(list=ls())

#install.packages("estimatr")
#install.packages("ri")
#install.packages("randomizr")
#install.packages("AER")
#install.packages("descrb")
#install.packages("nnet")
#install.packages("weights")

library(estimatr)
library(ri)
library(randomizr)
library(AER)
library(descr)
library(nnet)
library(weights)

options(scipen=999)

#Read in dataframe

dataset_complete_analysis<-readRDS("dataset_complete_analysis.rds")

#Table 1

Tab_1<-NA

Tab_1 <- data.frame(c( rep("Political interest", 1), rep("Knowledge male candidate",1), rep("Knowledge female candidate", 1), rep("Political ambition", 1)))

Tab_1$meanmale<- NA
Tab_1$sdmale <- NA
Tab_1$nmale <- NA
Tab_1$meanfemale <- NA
Tab_1$sdfemale <- NA
Tab_1$nfemale <- NA
Tab_1$dim <- NA
Tab_1$pvalue <- NA
colnames(Tab_1) <- c("Covariate", "Mean", "SD", "N", "Mean", "SD", "N", "Gender Gap", "P-Value")
dim(Tab_1)
Tab_1


#Political interest

dataset_complete_analysis$W1_pol_interest_mis[dataset_complete_analysis$pol_interest_mis==0]<-0
dataset_complete_analysis$W1_pol_interest_mis[dataset_complete_analysis$pol_interest_mis==1]<-NA

data_W1_pol_interest<-subset(dataset_complete_analysis,!is.na(W1_pol_interest_mis))

Tab_1[1,2]<-round(mean1<-mean(data_W1_pol_interest$pol_interest[data_W1_pol_interest$fem==0]), digits=2)

Tab_1[1,5]<-round(mean2<-mean(data_W1_pol_interest$pol_interest[data_W1_pol_interest$fem==1]), digits=2)

Tab_1[1,8]<-round(mean1-mean2, digits=2) 

tab1<-table(data_W1_pol_interest$fem, data_W1_pol_interest$pol_interest)

Tab_1[1,4]<-margin.table(tab1,1)[1]
Tab_1[1,7]<-margin.table(tab1,1)[2]

Tab_1[1,9]<-round(summary(lm_robust(data_W1_pol_interest$pol_interest~data_W1_pol_interest$fem))$coefficients[2,4] , digits=2)

Tab_1[1,3]<-round(sd1<-sd(data_W1_pol_interest$pol_interest[data_W1_pol_interest$fem==0]), digits=2)

Tab_1[1,6]<-round(sd2<-sd(data_W1_pol_interest$pol_interest[data_W1_pol_interest$fem==1]), digits=2)


#Knowledge male candidates

dataset_complete_analysis$W1_knowledge_male_mis<-1
dataset_complete_analysis$W1_knowledge_male_mis[dataset_complete_analysis$knowledge_male_mis==1]<-NA

data_W1_knowledge_male<-subset(dataset_complete_analysis,!is.na(W1_knowledge_male_mis))

Tab_1[2,2]<-round(mean1<-mean(data_W1_knowledge_male$W1_knowledge_male[data_W1_knowledge_male$fem==0]), digits=2)

Tab_1[2,5]<-round(mean2<-mean(data_W1_knowledge_male$W1_knowledge_male[data_W1_knowledge_male$fem==1]), digits=2)

Tab_1[2,8]<-round(mean1-mean2, digits=2) 

tab1<-table(data_W1_knowledge_male$fem, data_W1_knowledge_male$W1_knowledge_male)

Tab_1[2,4]<-margin.table(tab1,1)[1]
Tab_1[2,7]<-margin.table(tab1,1)[2]

Tab_1[2,9]<-round(summary(lm_robust(data_W1_knowledge_male$W1_knowledge_male~data_W1_knowledge_male$fem))$coefficients[2,4] , digits=2)

Tab_1[2,3]<-round(sd1<-sd(data_W1_knowledge_male$W1_knowledge_male[data_W1_knowledge_male$fem==0]), digits=2)

Tab_1[2,6]<-round(sd2<-sd(data_W1_knowledge_male$W1_knowledge_male[data_W1_knowledge_male$fem==1]), digits=2)


#Knowledge female candidates

dataset_complete_analysis$W1_knowledge_female_mis<-1
dataset_complete_analysis$W1_knowledge_female_mis[dataset_complete_analysis$knowledge_female_mis==1]<-NA

data_W1_knowledge_female<-subset(dataset_complete_analysis,!is.na(W1_knowledge_female_mis))

Tab_1[3,2]<-round(mean1<-mean(data_W1_knowledge_female$W1_knowledge_female[data_W1_knowledge_female$fem==0]), digits=2)

Tab_1[3,5]<-round(mean2<-mean(data_W1_knowledge_female$W1_knowledge_female[data_W1_knowledge_female$fem==1]), digits=2)

Tab_1[3,8]<-round(mean1-mean2, digits=2) 

tab1<-table(data_W1_knowledge_female$fem, data_W1_knowledge_female$W1_knowledge_female)

Tab_1[3,4]<-margin.table(tab1,1)[1]
Tab_1[3,7]<-margin.table(tab1,1)[2]

Tab_1[3,9]<-round(summary(lm_robust(data_W1_knowledge_female$W1_knowledge_female~data_W1_knowledge_female$fem))$coefficients[2,4] , digits=2)

Tab_1[3,3]<-round(sd1<-sd(data_W1_knowledge_female$W1_knowledge_female[data_W1_knowledge_female$fem==0]), digits=2)

Tab_1[3,6]<-round(sd2<-sd(data_W1_knowledge_female$W1_knowledge_female[data_W1_knowledge_female$fem==1]), digits=2)


#Running for office

dataset_complete_analysis$W1_running_for_office_mis<-1
dataset_complete_analysis$W1_running_for_office_mis[dataset_complete_analysis$running_for_office_mis==1]<-NA

data_W1_running_for_office<-subset(dataset_complete_analysis,!is.na(W1_running_for_office_mis))

Tab_1[4,2]<-round(mean1<-mean(data_W1_running_for_office$running_for_office[data_W1_running_for_office$fem==0]), digits=2)

Tab_1[4,5]<-round(mean2<-mean(data_W1_running_for_office$running_for_office[data_W1_running_for_office$fem==1]), digits=2)

Tab_1[4,8]<-round(mean1-mean2, digits=2) 

tab1<-table(data_W1_running_for_office$fem, data_W1_running_for_office$running_for_office)

Tab_1[4,4]<-margin.table(tab1,1)[1]
Tab_1[4,7]<-margin.table(tab1,1)[2]

Tab_1[4,9]<-round(summary(lm_robust(data_W1_running_for_office$running_for_office~data_W1_running_for_office$fem))$coefficients[2,4] , digits=2)

Tab_1[4,3]<-round(sd1<-sd(data_W1_running_for_office$running_for_office[data_W1_running_for_office$fem==0]), digits=2)

Tab_1[4,6]<-round(sd2<-sd(data_W1_running_for_office$running_for_office[data_W1_running_for_office$fem==1]), digits=2)

Tab_1


#########################################
# Analysis Field Experiment
#########################################

#Table 2

### Behavioral outcome ###


#Compliance

fieldexperiment <-subset(dataset_complete_analysis,!is.na(Treatment), select=c(Treatment, Block, Comply, Answer_invite, Application,age, married, Swiss_citizenship, year_uni, employment_father_1, employment_mother_1, edu_father_3,	edu_father_2, edu_father_1,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3, division_4,	division_6,	division_7, degree_1, degree_2, degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2, social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male))

set.seed(1234567)

table(dataset_complete_analysis$division_8)

Z <- fieldexperiment$Treatment

table(Z)

Block_var <- fieldexperiment$Block

D1 <- fieldexperiment$Comply

D2 <- fieldexperiment$Answer_invite

table(D1)

table(D2)

probs <- genprobexact(Z, blockvar=Block_var)
probs

w <- Z/probs + (1-Z)/(1-probs)

crosstab(Z, D1, prop.r = TRUE)

crosstab(Z, D2,  prop.r = TRUE)

crosstab(Z, D1, weight = w, prop.r = TRUE)

crosstab(Z, D2, weight = w, prop.r = TRUE)

fig_Comply<-round(crosstab(Z, D1, weight = w, prop.r = TRUE)$prop.row[2,2], digits=2)

fig_Answer_invite<-round(crosstab(Z, D2, weight = w, prop.r = TRUE)$prop.row[2,2], digits=2)


## ITT ##
table_ITT <-NA
table_ITT <- data.frame(c(rep("mean",1), rep("ITT", 2), rep("N", 1)), c(rep(c("No", "No", "Yes", ""))))
table_ITT$ITT_behav <- NA
table_ITT$cii_behav <- NA
table_ITT$cis_behav <- NA
table_ITT$ITT_att <- NA
table_ITT$cii_att <- NA
table_ITT$cis_att <- NA
colnames(table_ITT) <- c("Estimate", "Covariate_adjustment", "ITT", "cil", "ciu", "ITT", "cil", "ciu")
dim(table_ITT)
table_ITT

set.seed(1234567)
Z <- fieldexperiment$Treatment
Blockvar<-fieldexperiment$Block
Y <- fieldexperiment$Application
X <- as.matrix(fieldexperiment[,6:47])
D1<-fieldexperiment$Comply
D2<-fieldexperiment$Answer_invite

table(D1)
table(D2)
table(Y)
probs <- genprobexact(Z, blockvar=Blockvar)
probs

w <- Z/probs + (1-Z)/(1-probs)

library(descr)

table_ITT[1,3]<-round(crosstab(Z, Y, weight = w, prop.r = TRUE)$prop.row[1,2], digits=2)


probs <- genprobexact(Z, blockvar=Blockvar)
probs
ate <- estate(Y,Z,prob=probs)
round(ate, 2) 
table_ITT[2,3] <- round(ate, digits=2)

set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)
Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,prob=probs)
dispdist(distout,ate) 
names(dispdist(distout,ate))
round(dispdist(distout,ate)$sd, 2) 

Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,prob=probs)
ate
dispdist(distout,ate)
table_ITT[2,4:5] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)


# Covariate-adjustment

probs <- genprobexact(Z, blockvar=Blockvar)
probs
ate <- estate(Y,Z,X,prob=probs)
round(ate, 2) 
table_ITT[3,3] <- round(ate, digits=2)

set.seed(1234567)

perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)

Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,X,prob=probs)
dispdist(distout,ate) 
round(dispdist(distout,ate)$sd, 2) 

Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,X,prob=probs)
ate
dispdist(distout,ate) 
round(dispdist(distout,ate)$quantile, 2)
table_ITT[3,4:5] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)

table_ITT[4,3] <- table(fieldexperiment$Treatment)[1]+table(fieldexperiment$Treatment)[2]


table_ITT


### Attitudinal outcome ("could see myself running for office") ###

fieldexperiment1 <-subset(dataset_complete_analysis,!is.na(Treatment), select=c(Treatment, Block, Comply, Answer_invite, running_for_office_w2, age, married,  Swiss_citizenship, year_uni, employment_father_1, employment_father_3, employment_mother_1, employment_mother_3, edu_father_4,	edu_father_3,	edu_father_2, edu_father_1,	edu_mother_4,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3,	division_4,	division_5,	division_6,	division_7, degree_1, degree_2,	degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2,	social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male, knowledge_male_mis, knowledge_female_mis))

fieldexperiment1 <-subset(fieldexperiment1,!is.na(running_for_office_w2), select=c(Treatment, Block, Comply, Answer_invite, running_for_office_w2, age, married, Swiss_citizenship, year_uni, employment_father_1, employment_mother_1, edu_father_3,	edu_father_2, edu_father_1,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3, division_4,	division_6,	division_7, degree_1, degree_2, degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2, social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male))

table(dataset_complete_analysis$running_for_office)

Z <- fieldexperiment1$Treatment
Blockvar <- fieldexperiment1$Block
Y <- fieldexperiment1$running_for_office_w2
X <- as.matrix(fieldexperiment1[,6:47])
D1 <- fieldexperiment1$Comply
D2 <- fieldexperiment1$Answer_invite

table(Z)
table(D1)

probs <- genprobexact(Z, blockvar=Blockvar)
probs

w <- Z/probs + (1-Z)/(1-probs)

table_ITT[1,6]<-round(crosstab(Z, Y, weight = w, prop.r = TRUE)$prop.row[1,2], digits=2)


ate <- estate(Y,Z,prob=probs)
round(ate, 2) 
table_ITT[2,6] <- round(ate, digits=2)

set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)
Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,prob=probs)
dispdist(distout,ate) 
round(dispdist(distout,ate)$sd, 2) 


Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,prob=probs)
ate
dispdist(distout,ate) 
round(dispdist(distout,ate)$quantile, 2)
table_ITT[2,7:8] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)


# Covariate-adjustment

probs <- genprobexact(Z, blockvar=Blockvar)
probs
ate <- estate(Y,Z,X,prob=probs)
ate 
round(ate, 2)
table_ITT[3,6] <- round(ate, digits=2)


set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)

Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,X,prob=probs)
dispdist(distout,ate) 

Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,X,prob=probs)
ate
dispdist(distout,ate) 
round(dispdist(distout,ate)$quantile, 2)
table_ITT[3,7:8] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)

table_ITT[4,6] <- table(fieldexperiment1$Treatment)[1]+table(fieldexperiment1$Treatment)[2]


table_ITT



#Table 3


## CACE ##

table_CACE <-NA
table_CACE <- data.frame(c(rep("mean",1), rep("CACE", 2), rep("N", 1)), c(rep(c("No", "No", "Yes", ""))))
table_CACE$CACE_behav <- NA
table_CACE$cil_behav <- NA
table_CACE$cih_behav <- NA
table_CACE$CACE_att <- NA
table_CACE$cil_att <- NA
table_CACE$cih_att <- NA
colnames(table_CACE) <- c("Estimate", "Covariate_adjustment", "CACE", "cil", "ciu", "CACE", "cil", "ciu")
dim(table_CACE)
table_CACE


fieldexperiment <-subset(dataset_complete_analysis,!is.na(Treatment), select=c(Treatment, Block, Comply, Answer_invite, Application,age, married, Swiss_citizenship, year_uni, employment_father_1, employment_mother_1, edu_father_3,	edu_father_2, edu_father_1,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3, division_4,	division_6,	division_7, degree_1, degree_2, degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2, social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male))

set.seed(1234567)
Z <- fieldexperiment$Treatment
Blockvar<-fieldexperiment$Block
Y <- fieldexperiment$Application
X <- as.matrix(fieldexperiment[,6:47])
D1<-fieldexperiment$Comply
D2<-fieldexperiment$Answer_invite


table(D1)
table(D2)
table(Y)
probs <- genprobexact(Z, blockvar=Blockvar)
probs

w <- Z/probs + (1-Z)/(1-probs)


table_CACE[1,3]<-round(crosstab(Z, D1, weight = w, prop.r = TRUE)$prop.row[2,2], digits=2)

table_CACE[1,6]<-round(crosstab(Z, D1, weight = w, prop.r = TRUE)$prop.row[2,2], digits=2)

CACE1 <- ivreg(Y~D1 | Z, weight=w)
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_CACE[2,3] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_CACE[2,4] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_CACE[2,5] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)


CACE2 <- ivreg(Y~D1 + X| Z + X, weight=w)
summary(CACE2)
coeftest(CACE2,vcovHC(CACE2, type="HC2"))
round(coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], 2)
round(coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], 2)
table_CACE[3,3] <- round(coeftest(CACE2,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_CACE[3,4] <- round(coeftest(CACE2,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], digits=2)
table_CACE[3,5] <- round(coeftest(CACE2,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], digits=2)

table_CACE[4,3] <- table(fieldexperiment$Treatment)[1]+table(fieldexperiment$Treatment)[2]


## CACE Att ##

fieldexperiment1 <-subset(dataset_complete_analysis,!is.na(Treatment), select=c(Treatment, Block, Comply, Answer_invite, running_for_office_w2, age, married,  Swiss_citizenship, year_uni, employment_father_1, employment_father_3, employment_mother_1, employment_mother_3, edu_father_4,	edu_father_3,	edu_father_2, edu_father_1,	edu_mother_4,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3,	division_4,	division_5,	division_6,	division_7, degree_1, degree_2,	degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2,	social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male, knowledge_male_mis, knowledge_female_mis))

fieldexperiment1 <-subset(fieldexperiment1,!is.na(running_for_office_w2), select=c(Treatment, Block, Comply, Answer_invite, running_for_office_w2, age, married, Swiss_citizenship, year_uni, employment_father_1, employment_mother_1, edu_father_3,	edu_father_2, edu_father_1,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3, division_4,	division_6,	division_7, degree_1, degree_2, degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2, social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male))


Z <- fieldexperiment1$Treatment
Blockvar <- fieldexperiment1$Block
Y <- fieldexperiment1$running_for_office_w2
X <- as.matrix(fieldexperiment1[,6:47])
D1 <- fieldexperiment1$Comply
D2 <- fieldexperiment1$Answer_invite

table(Z)
table(D1)

probs <- genprobexact(Z, blockvar=Blockvar)
probs

w <- Z/probs + (1-Z)/(1-probs)


CACE1 <- ivreg(Y~D1 | Z, weight=w)
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_CACE[2,6] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_CACE[2,7] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_CACE[2,8] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)


CACE2 <- ivreg(Y~D1 + X| Z + X, weight=w)
summary(CACE2)
coeftest(CACE2,vcovHC(CACE2, type="HC2"))
round(coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], 2)
round(coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], 2)
table_CACE[3,6] <- round(coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,1], digits=2)
table_CACE[3,7] <- round(coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], digits=2)
table_CACE[3,8] <- round(coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE2,vcovHC(CACE2, type="HC2"))[2,2], digits=2)

table_CACE[4,6] <- table(fieldexperiment1$Treatment)[1]+table(fieldexperiment1$Treatment)[2]

#Table 3
table_CACE


#Table 4
#media

table_media <-NA
table_media <- data.frame(c(rep("mean", 1), rep("comply",1), rep("ITT", 1), rep("ITT", 1), rep("CACE", 1), rep("CACE", 1), rep("N", 1)), c("", "", "No", "Yes", "No", "Yes", ""))
table_media$ITT_media <- NA
table_media$cil_media <- NA
table_media$cih_media <- NA
table_media$ITT_online <- NA
table_media$cil_online <- NA
table_media$cih_online <- NA
table_media$ITT_knowledge <- NA
table_media$cil_knowledge <- NA
table_media$cih_knowledge <- NA
colnames(table_media) <- c("Estimate", "Covariate_adjustment", "ITT", "cil", "ciu", "ITT", "cil", "ciu", "ITT", "cil", "ciu")
table_media


fieldexperiment3 <- subset(dataset_complete_analysis,!is.na(Treatment), select=c(Treatment, Block, Comply, Answer_invite, campaign_media_w2, age, married,  Swiss_citizenship, year_uni, employment_father_1,  employment_mother_1, edu_father_4,  edu_father_3,	edu_father_2,	edu_father_1,	edu_mother_4,	edu_mother_3,	edu_mother_2,	edu_mother_1, division_1, division_2,	division_3,	division_4,	division_5,	division_6,	division_7,	degree_1, degree_2,	degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2,	social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male, knowledge_male_mis, knowledge_female_mis))
fieldexperiment3[21]
fieldexperiment3 <- subset(fieldexperiment3,!is.na(campaign_media_w2), select=c(Treatment, Block, Comply, Answer_invite, campaign_media_w2, age, married, Swiss_citizenship, year_uni, employment_father_1, employment_mother_1, edu_father_3,	edu_father_2, edu_father_1,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3, division_4,	division_6,	division_7, degree_1, degree_2, degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2, social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male))


Z <- fieldexperiment3$Treatment
Blockvar <- fieldexperiment3$Block
Y <- fieldexperiment3$campaign_media_w2
X <- as.matrix(fieldexperiment3[,6:47])
D1 <- fieldexperiment3$Comply
D2 <- fieldexperiment3$Answer_invite

probs <- genprobexact(Z, blockvar=Blockvar)
probs

w <- Z/probs + (1-Z)/(1-probs)

table_media[1,3]<-round(weighted.mean(Y, weight = w[Z=0]), digits=2)

table_media[2,3]<-round(crosstab(Z, D1, weight = w, prop.r = TRUE)$prop.row[2,2], digits=2)


ate <- estate(Y,Z,prob=probs)
round(ate, 2) 
table_media[3,3] <- round(ate, digits=2)

set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)
Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,prob=probs)
dispdist(distout,ate) 
round(dispdist(distout,ate)$sd, 2)


Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,prob=probs)
ate
dispdist(distout,ate)
round(dispdist(distout,ate)$quantile, 2)
table_media[3,4:5] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)


# Covariate-adjustment

probs <- genprobexact(Z, blockvar=Blockvar)
probs
ate <- estate(Y,Z,X,prob=probs)
ate 
round(ate, 2) 
table_media[4,3] <- round(ate, digits=2)


set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)

Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,X,prob=probs)
dispdist(distout,ate) 

Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,X,prob=probs)
ate
dispdist(distout,ate) 
round(dispdist(distout,ate)$quantile, 2)
table_media[4,4:5] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)

# CACE

w <- Z/probs + (1-Z)/(1-probs)
CACE1 <- ivreg(Y~D1 | Z, weight=w)
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_media[5,3] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_media[5,4] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_media[5,5] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)


CACE1 <- ivreg(Y~D1 + X| Z + X, weight=w)
summary(CACE1)
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_media[6,3] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_media[6,4] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_media[6,5] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)

table_media[7,3] <- round(table(fieldexperiment3$Treatment)[1]+table(fieldexperiment3$Treatment)[2], digits=0)


#follow campaign online

table(dataset_complete_analysis$campaign_online_w2)

fieldexperiment4 <-subset(dataset_complete_analysis,!is.na(Treatment), select=c(Treatment, Block, Comply, Answer_invite, campaign_online_w2, age, married,  Swiss_citizenship, year_uni, employment_father_1,  employment_mother_1, edu_father_4,  edu_father_3,  edu_father_2,	edu_father_1,	edu_mother_4,	edu_mother_3,	edu_mother_2,	edu_mother_1, division_1, division_2,	division_3,	division_4,	division_5,	division_6,	division_7,	degree_1, degree_2,	degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2,	social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male, knowledge_male_mis, knowledge_female_mis))
fieldexperiment4[21]
fieldexperiment4 <-subset(fieldexperiment4,!is.na(campaign_online_w2), select=c(Treatment, Block, Comply, Answer_invite, campaign_online_w2, age, married, Swiss_citizenship, year_uni, employment_father_1, employment_mother_1, edu_father_3,	edu_father_2, edu_father_1,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3, division_4,	division_6,	division_7, degree_1, degree_2, future_per,	career_goals_1, career_entry_1, children_1, children_2, social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male))
dim(fieldexperiment4)

Z <- fieldexperiment4$Treatment
Blockvar <- fieldexperiment4$Block
Y <- fieldexperiment4$campaign_online_w2
X <- as.matrix(fieldexperiment4[,6:45])
D1 <- fieldexperiment4$Comply
D2 <- fieldexperiment4$Answer_invite

table(Y)
table(Z)
table(Z, D1)

probs <- genprobexact(Z, blockvar=Blockvar)
probs

w <- Z/probs + (1-Z)/(1-probs)

table_media[1,6]<-round(weighted.mean(Y, weight = w[Z=0]), digits=2)

table_media[2,6]<-round(crosstab(Z, D1, weight = w, prop.r = TRUE)$prop.row[2,2], digits=2)


ate <- estate(Y,Z,prob=probs)
round(ate, 2) 
table_media[3,6] <- round(ate, digits=2)

set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)
Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,prob=probs)
dispdist(distout,ate) 
round(dispdist(distout,ate)$sd, 2)


Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,prob=probs)
ate
dispdist(distout,ate) 
round(dispdist(distout,ate)$quantile, 2)
table_media[3,7:8] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)


# Covariate-adjustment

probs <- genprobexact(Z, blockvar=Blockvar)
probs
ate <- estate(Y,Z,X,prob=probs)
ate 
round(ate, 2) 
table_media[4,6] <- round(ate, digits=2)


set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)

Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,X,prob=probs)
dispdist(distout,ate)

Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,X,prob=probs)
ate
dispdist(distout,ate)
round(dispdist(distout,ate)$quantile, 2)
table_media[4,7:8] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)


# CACE

w <- Z/probs + (1-Z)/(1-probs)
CACE1 <- ivreg(Y~D1 | Z, weight=w)
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_media[5,6] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_media[5,7] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_media[5,8] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)


CACE1 <- ivreg(Y~D1 + X| Z + X, weight=w)
summary(CACE1)
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_media[6,6] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_media[6,7] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_media[6,8] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)

table_media[7,6] <- round(table(fieldexperiment4$Treatment)[1]+table(fieldexperiment4$Treatment)[2], digits=0)


#Knowledge

fieldexperiment_knowledge<-subset(dataset_complete_analysis,!is.na(Treatment), select=c(Treatment, Block, Comply, Answer_invite, W2_knowledge, age, married,  Swiss_citizenship, year_uni, employment_father_1,  employment_mother_1, edu_father_4,  edu_father_3,  edu_father_2,  edu_father_1,	edu_mother_4,	edu_mother_3,	edu_mother_2,	edu_mother_1, division_1, division_2,	division_3,	division_4,	division_5,	division_6,	division_7,	degree_1, degree_2, degree_3,	future_per,	career_goals_1, career_entry_1, children_1, children_2,	social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male, knowledge_male_mis, knowledge_female_mis))

fieldexperiment_knowledge <-subset(fieldexperiment_knowledge,!is.na(W2_knowledge), select=c(Treatment, Block, Comply, Answer_invite, W2_knowledge, age, married, Swiss_citizenship, year_uni, employment_father_1, employment_mother_1, edu_father_3,	edu_father_2, edu_father_1,	edu_mother_3,	edu_mother_2, edu_mother_1, division_1, division_2,	division_3, division_4,	division_6,	division_7, degree_1, degree_2, future_per,	career_goals_1, career_entry_1, children_1, children_2, social_skills,	presentation_skills,	networking_skills,	running_for_office,	politicians	,	gender_childcare,	gender_career,	gender_election	, gender_work,	gender_politics, gender_children,	econ_interest,	pol_interest,	pol_knowledge_self_report, W1_knowledge_female, W1_knowledge_male))
table(fieldexperiment_knowledge)
Z <- fieldexperiment_knowledge$Treatment
Blockvar <- fieldexperiment_knowledge$Block
Y <- fieldexperiment_knowledge$W2_knowledge
X <- as.matrix(fieldexperiment_knowledge[,6:36])
D1 <- fieldexperiment_knowledge$Comply
D2 <- fieldexperiment_knowledge$Answer_invite

dim(fieldexperiment_knowledge)

table(Y)
table(Z)

probs <- genprobexact(Z, blockvar=Blockvar)
probs

w <- Z/probs + (1-Z)/(1-probs)

table_media[1,9]<-round(weighted.mean(Y, weight = w[Z=0]), digits=2)

table_media[2,9]<-round(crosstab(Z, D1, weight = w, prop.r = TRUE)$prop.row[2,2], digits=2)


ate <- estate(Y,Z,prob=probs)
round(ate, 2) 
table_media[3,9] <- round(ate, digits=2)

set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)
Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,prob=probs)
dispdist(distout,ate) 
round(dispdist(distout,ate)$sd, 2) 


Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,prob=probs)
ate
dispdist(distout,ate) 
round(dispdist(distout,ate)$quantile, 2)
table_media[3,10:11] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)


# Covariate-adjustment

probs <- genprobexact(Z, blockvar=Blockvar)
probs
ate <- estate(Y,Z,X,prob=probs)
ate 
round(ate, 2) 
table_media[4,9] <- round(ate, digits=2)


set.seed(1234567)
perms <- genperms(Z,maxiter=5000, blockvar=Blockvar)

Ys <- genouts(Y,Z,ate=0)
distout <- gendist(Ys,perms,X,prob=probs)
dispdist(distout,ate) 

Ys <- genouts(Y,Z,ate=ate)
distout <- gendist(Ys,perms,X,prob=probs)
ate
dispdist(distout,ate) 
round(dispdist(distout,ate)$quantile, 2)
table_media[4,10:11] <- round(dispdist(distout,ate, quantiles = c(0.025, 0.975))$quantile, digits=2)



# CACE

w <- Z/probs + (1-Z)/(1-probs)
CACE1 <- ivreg(Y~D1 | Z, weight=w)
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_media[5,9] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_media[5,10] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_media[5,11] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)


CACE1 <- ivreg(Y~D1 + X| Z + X, weight=w)
summary(CACE1) 
coeftest(CACE1,vcovHC(CACE1, type="HC2"))
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], 2)
table_media[6,9] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1], digits=2)
table_media[6,10] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] - qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)
table_media[6,11] <- round(coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,1] + qnorm(0.975)*coeftest(CACE1,vcovHC(CACE1, type="HC2"))[2,2], digits=2)


table_media[7,9] <- round(table(fieldexperiment_knowledge$Treatment)[1]+table(fieldexperiment_knowledge$Treatment)[2], digits=0)

#Table 4
table_media


#Compile log file

sink("log_results.txt")

#Table 1
print(Tab_1)

#Table 2
print(table_ITT)

#Table 3
print(table_CACE)

#Table 4
print(table_media)

sink()
